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Starting from the von Neumann-Maxweh equations for the Wigner quasi-probabihty distribution and for the 
self-consistent electric field, the quantum analog of the classical single-wave model has been derived. The 
linear stability of the quantum single-wave model has been studied, and periodic in time patterns have been 
found both analytically and numerically. In addition, some features of quantum chaos have been detected in 
the unstable region in parameter space. Further, a class of standing-wave solutions of the quantum single- 
wave model has also been found, which have been observed to behave as stable solitary-wave structures. 
The analytical results have been finally compared to the exact system dynamics obtained by solving the 
corresponding equations in Schrodinger representation numerically. 
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I. INTRODUCTION 

Conventional plasma physics describes dynamical and 
equilibrium processes in systems characterized by high 
temperatures and relatively low densities, at which quan- 
tum effects have practically no influence. Recent techno- 
logical progress on miniaturized semiconductor devices 
and nanoscale structures has stimulated the interest in 
potential applications of plasma physics, where the quan- 
tum nature of the constituting objects plays an essential 
role. 

Both fusion and space plasmas are characterized by 
states of high temperatures and low densities, for which 
quantum effects are negligible and are usually not taken 
into account. The classical example, where plasma and 
quantum effects mutually coexist is the electron gas in 
metals at room temperature. The electron population 
in ordinary metals is globally neutralized by the lattice 
ions and therefore can be regarded as a classical exam- 
ple of a plasma system. However, its basic properties 
are mainly determined by the presence of the regular ion 
lattice, so that typical plasma effects are only correction 
of higher order. In recent years, substantial progress has 
been made in preparation and manipulation of metal- 
lic nanostructures consisting of a small number of atoms 
(typically 10 — 10^)^. Since an underlying ionic lattice 
for such structures does not exist, the electron dynamics 
is dominated by plasma effects, at least for systems with 
sufficiently large spatial dimensions. This makes metallic 
nanostructures the ideal object to study equilibrium and 
evolution properties of quantum plasmas. 
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Additional important applications of quantum plasmas 
arise from semiconductor physics and astrophysics^Ti^. 
The basic factor dictating the necessity of utilizing quan- 
tum mechanical description in semiconductors is the 
great degree of miniaturization of today's electronic com- 
ponents to the extent that the de Broglie wavelength of 
the charge carriers can be comparable to the spatial vari- 
ation of the doping profiles. In some astrophysical ob- 
jects under extreme conditions of temperature and den- 
sity, such as white dwarf stars, the density can become 
some ten order of magnitudes larger than that of ordinary 
solids. Because of such large densities, a white dwarf can 
be as hot as a fusion plasma (10* K), but still behave as 
a quantum-mechanical object. 

When quantum effects start playing a significant role, 
the well-known picture from conventional plasma physics 
becomes even more complicated, because an additional 
length scale has to be introduced. This is the de Broglie 
wavelength of a charged particle, Xb = h/[mvT), where 
vt is its thermal velocity. The latter represents approxi- 
mately the spatial extension of the particle wave function 
and obviously the larger it is, the more significant quan- 
tum effects are. It is also clear that quantum behavior 
will be reached more easily for the electrons than for the 
ions, due to the large mass difference. Indeed, in all prac- 
tical situations, even in the most extreme ones, the ion 
dynamics remains always classical, while only the elec- 
trons require quantum-mechanical description. 

Another important area uncovering broad vistas for 
applications of the quantum plasma dynamics approach 
is the free electron laser (FEL). If the photon momentum 
recoil is greater than, or of the order of the momentum 
spread in the electron beam traversing the undulator, 
quantum effects contribute substantially and cannot be 
neglected''. 

Due to the tremendous number of potential applica- 
tions ranging from recent high technology developments 
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to high energy physics and astrophysics, there has been a 
growing interest in studying new aspects of dense quan- 
tum plasmas^. 

The purpose of the present paper is to extend the clas- 
sical single-wave model^ii^ to the quantum case. Follow- 
ing del-Castillo-Negrete^ we give in the next Section a 
brief sketch of the derivation of the quantum single- wave 
model. In Sections IIIII and IIVI the perturbation solution 
of the quantum single-wave model in Schrodinger pic- 
ture for a mixed state with different energies has been 
constructed, and in Section |V] the amplitude equation 
describing the nonlinear wave interaction between eigen- 
modes has been obtained. Further, in Section IVll a class 
of solutions of the quantum single- wave model in the form 
of nonlinear standing waves has been found, much resem- 
bling solitary wave structures. In Section IVIII the com- 
parison between the theory and direct numerical simu- 
lation of the underlying equations in Schrodinger repre- 
sentation has been presented. Finally, in Section IVlIII we 
draw some conclusions. 



II. SKETCH OF DERIVATION OF THE QUANTUM 
SINGLE-WAVE MODEL 



The von Neumann-Maxwell system of equations is cast 
now in the form 

OtF + PdxF 



dl^ = J dPF{X,P;T), 



(Fo + F), 
(5) 

(6) 



where the new dimensionless Planck's constant is given 
by the expression 



flWp 

hi- 



mvf. 



(7) 



Using the formal Fourier representation of the Wigner 
quasi-probability distribution in momentum space, and 
performing simple algebraic manipulations, it is possible 
to cast the von Neumann equation in an equivalent form 

OtF + PdxF 



We start with the von Neumann-Maxwell equations in 
one spatial dimension 

dt,F+^dx,F 
m 



te 



ih , 



ih „ 

X, - y5p, 



{Fo 



J dPiF(Xi,Fi 



eiVo 
F{Xi,Pi;h) 



h). 



(1) 
(2) 



is the Wigner quasi- 
m and e are the elec- 



Here, Fo(Pi) 

probability distribution functior 
tron rest mass and charge, respectively, and eo is the 
vacuum permittivity. Furthermore, Nq is the density of 
the background ions, while all independent variables with 
index "1" refer to the real time and phase-space vari- 
ables. The latter is adopted only for the sake of clarity 
in distinguishing the final notation for the independent 
variables in a sequence of scaling transformations per- 
formed throughout the derivation. Note that the quasi- 
probability distribution function F{Xi, Pi;ti) is the de- 
parture from the equilibrium one Fo(Pi). It is convenient 
to nondimensionalize the variables in the above Eqs. ([T]) 
and ^ according to 



T 



X 



P 



(3) 



2TTh 



- J dAdPie'^(^-^i) 



$( X - ^) - <i>(^X + [Fo(Pi) + F{X, P,-T)], 

. 

which sometimes proves to be more convenient for direct 
analysis. 

The linear stability of the equilibrium electron distri- 
bution Fq{P) is now determined by the roots of the dis- 
persion function 

X>(fc,c) = 



1 f dP 



fiifc3 J P- 



F^[ P 



hik 



FoiP 



hik 



(9) 

where k is the wave number, and c — uj/k. The quantum 
correction to the dispersion function (|9]) as compared to 
the classical limiting case, when hi — 0, is of second 
order in hi 



V{k,c)^V,iik,c)~^ J 



hi f dP d^FoiP) 



P~c dP3 



(10) 



where Vc is a velocity characteristic for the system, which 
is to a large extent arbitrary, while ujp is the plasma fre- 
quency 



(4) 



where 



2),,(fc,c) = l-l J 



1 f dP dFo(F) 



P-c dP 



fll) 
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is the classical dispersion function. Therefore, for sufB- 
ciently small hi one should not expect significant modi- 
fications in the linear stability properties as described by 
the roots of the classical dispersion function ([TT|) . 

Following del-Castillo-Negrete^, we assume the exis- 
tence of a marginally stable state, prescribed by a station- 
ary inflection point of the equilibrium distribution Fq (P) . 
Since the analysis does not add anything new apart from 
taking into account the "quantum" Vlasov equation ^ 
[or its equivalent form ([5])] instead of the classical one, 
we shall omit details of the calculation here. The final 
result is 



h 



(p{x;t) ^a{t)e''' +a*{t)e-", 
where the operator averaging is denoted by 

oo 27r 
(...) = ^/dp/d..... 




(14) 



(15) 



The dimensionless time t, coordinate x and momen- 
tum p, as well as the single wave amplitude a and the 
quasi-probability distribution function / are appropri- 
ately nondimensionalized^ . In addition, the new dimen- 
sionless Planck's constant H is defined as 



^_ 2^^l7l 
mVcL\X\ ' 

where L is the spatial length of the system, and 



(16) 



(17) 



In the above definition p7)) . (A:*,c*) characterizes the 
marginally stable state in the (fc, c)-space, and A is a 
parameter measuring the difference between the domain 
length L and the wavelength 27r/fc, of the inflection point 
mode. To simplify notations in what follows, we shall 
drop the tilde over the h as defined by Eq. p^ . 

The linear stability analysis of the quantum single- 
wave model as described by the von-Neumann equation 
(IT^ , coupled with the equations for the single field mode 
(IT^ and is presented in Appendix 



III. THE SCHRODINGER PICTURE OF THE 
QUANTUM SINGLE-WAVE MODEL 

The von Neumann equation (also known as Liouville- 
von Neumann equation) describes how the Wigner 
quasi-probability distribution evolves in time, while the 



Schrodinger equation describes how pure states with sup- 
port in configuration space evolve in time. In fact, the 
two equations are equivalent, in the sense that either can 
be derived from the other^^. More convenient for the sub- 
sequent analysis is the Schrodinger form of the system of 
equations (IT^ - (HH), which can be written as 



(/3(x;t) =a(t)e"-f a*(i)e" 



da 

a— — h ila 
at 



(18) 



(19) 



(20) 



where now the averaging specified by Eq. (jlSp involves 
integration on the spatial variable x only. It is worthwhile 
mentioning that an alternative quantum hydrodynamics 
approach (QHD) is possible^^'-'^*, however for the present 
purposes we prefer to confine ourselves to the Schrodinger 
picture. 

Taking into account the intrinsic periodicity of the sys- 
tem in the spatial variable a;, we represent the wave func- 
tion in the form 



^{x-t) 



OO 

E 



The dynamics of the Fourier harmonics 
by the following equation 



(21) 



is governed 



ihBt 



h m 



m+l, 



(22) 



which should be supplemented by the equation 

oo 

{adt+il)a^i (23) 

n— — oo 

describing the evolution of the amplitude of the self- 
consistent electric field single-wave mode. 

Let ko be a given harmonic number. Since direct inter- 
action between waves as described by Eq. (^5]) involves 
only closest neighbors, it is straightforward to verify that 
an exact solution of the system of equations and 
is of the form 



. .fco-2, fco,fco + 2, 



(24) 

and '^l? = otherwise. In addition, the amplitude of 
the single wave mode vanishes (a'^'^' = O) as should be 
expected. This solution represents a superposition of non 
interacting plane waves, which remain uncoupled due to 
the fact that the nonzero wave harmonics are of either 
even or odd wave number, respectively. 

To begin with, we consider the case, where the only 
wave with nonzero amplitude is the one with a prescribed 
harmonic number s {As 7^ 0), while all others vanish. 
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The sohition we shall be seeking can be expressed as fol- 
lows 



where 



^(0) _ Oms -ihm^t/2+iip„ 

™ ATT 
V ^TT 



(25) 
(26) 



and e is a formal small parameter. Thus, we obtain 



ihdt 



A„i = -e(Q;A„i_i + a*A„i+i), (27) 



for m 7^ s — 1 or m 7^ s + 1, and 



ihdt + -{s-if 



A:_i=a*(0)*+,(«*A:_2 + aA:), 

(28) 

A^+i = -a^f^ - e{aKs + Q*Ks+2), 

(29) 

{adt + il)a = z(*i°)A:_i + *i°)*A,+i) +ie^ A„A,*_i, 

n 

(30) 

It is convenient to introduce new variables defined ac- 
cording to the expressions 

i = *1°'a:_i + *f *A,+i, c - *i°^A:_i - *i°)*A,+i, 

(31) 

where 

A,_i=^vi/(o)(^*+C*), A,+i=^vi/(o)(e-C), (32) 

Simple algebraic manipulations of Eq. (f28| - (pO| yield 
the following equations for the new variables ^ and C, 

{ihdt-h^s)^+^C^^W^. (33) 
{ihdt-h^s)C+^^ = ^+eW^, (34) 

where 

= a**(°)A*_2 + a*i°^A: - a^f>Ks - a*^f>Ks+2, 

(35) 

= a**i°' A:_2 + a: + a*(°)*A, -f a**(")*A,+2. 

(36) 

Eliminating the auxiliary variable C, we obtain a single 
equation for ^ 



(aSt -I- iO(*9t - hsf - —{adt + ^ 

4 ZTT 



-£^A„A,:_i+|(^9t + *0 



(i9t - hs)W^ - -Wc 



(37) 

which should be supplemented by Eq. pop . written as 
(aa* +zZ)a = zC + ieXl^»^"-i- (3^) 



IV. PERTURBATIVE SOLUTION OF THE EVOLUTION 
EQUATIONS 

Following the standard perturbation approximation 
procedurei^, we represent all dynamical variables as a 
series expansion in the formal small parameter e. In or- 
der to eliminate resonant (secular) terms appearing in 
successive orders, we define different time scales as pre- 
scribed by the method of multiple scales^^. To zero order, 
the perturbation equations (1371) and (I38p can be written 
as 



ta ■ 
{adt + il){idt - hsf - —{adt +il) + 7r 

4 /TT 



Co =0, 
(39) 

{adt + il)aQ=i^Q. (40) 
The general solution of Eq. can be expressed as 

3 

Co(t,ii,i2,...) = E^/5(^i'^2,...)e'"''*, (41) 

p=l 

where ti = et, t2 = e^t, and so on are slow times, and 
are the three roots of the characteristic equation 

{au + l){u + hsf {(Tuj + l) + — =0. (42) 

4 27r 

The zero-order amplitude of the single- wave mode ao can 
be determined immediately from Eq. (PU]) to be 



(43) 



From Eq. ([33]) for the auxiliary variable ^Oi we obtain 

3 



(44) 



/3=1 



For the perturbed harmonics of the wave function we can 
now write 

Ai"_\ - E A^-^Ble-<^'^+'''"/^y+''^^ , (45) 



where the coefficients in the above expressions are given 
by 

A^"' = 1 + + ^s), A^+^ = 1 - ^{ujp + hs). 

(47) 

Note that in zero order all other harmonics as described 
by Eq. dUl) vanish (^A^^ = 0, m^s -1,3 + 1^. 
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Even in the lowest order an interesting feature of quan- 
tum single-wave model can be observed. Two adjacent 
harmonics (with wave numbers s — 1 and s -I- 1) are gen- 
erated on either side of the initial isolated harmonic with 
wave number s. This excitation is due to the closest 
neighbor interaction, characterized by the evolution law 
of the amplitude of the self-consistent electric potential. 
As it will become clear from the subsequent exposition, 
once harmonics s — 1 and s + 1 are generated, they give 
rise to their closest neighbors s — 2 and s + 2, respectively, 
and so on until all the spectrum becomes eventually pop- 
ulated. This phenomenon has also been detected in the 
numerical simulations. 

Let us proceed with the first order in the formal expan- 
sion parameter e. Since the quantity satisfies equation 
similar to Eq. (|39l). without loss of generality it can be 
assumed = 0. This immediately leads to 



Ci = 0, 



ai = 0. 



(48) 



Thus, the first order corrections to the (s — l)-st and 
(s + l)-st harmonics of the wave function also vanish 



A 



(1) 



0, 



A 



(1) 

s+l 



0. 



(49) 



Additional consequence of the vanishing first-order so- 
lution ^1 = is the fact that the amplitudes Bp intro- 
duced in Eq. (|4T|) do not depend on the time scale ti 
(dtAo = 0). 

Taking into account Eq. ((27)) and the fact that the 
arbitrary phase tps may depend on slower time scales, we 
obtain 



ihdt 



a*A^''^ 



(50) 

The right-hand-side of the above equation contains reso- 
nant (secular) terms, which can be eliminated by equat- 
ing them to the second term on the left-hand-side. Thus, 
we obtain an equation describing the phase dynamics of 
the basic harmonic '3/1°^ [see Eq. In the case, where 

all three roots of the characteristic equation (|42p are real, 
the amplitude dependent phase shift is governed by the 
following equation 



dt^s 



2vry \Bp\' _ 



(51) 



amplitude of the basic harmonic, we obtain 



Ad) = i /^^' ^^^^ 



A. 



(-) 



A 



(+) 



(53) 



X e 



The prime over the summation symbol in Eq. (j53p im- 
plies that terms for which o;^ = w* are excluded. Note 
that such terms have already been taken into account 
in writing the amplitude-phase equations (j5ip and (j52p . 
respectively. 

As already mentioned above, a new harmonic pair with 
wave numbers s — 2 and s -I- 2 is generated in first order. 
From Eq. (P7)) we obtain 



ihdt 



A 



(1) 

.s-2 



* A (0) 



(54) 



A^ll, = -aoAfl. (55) 



The solutions of the last two equations arc regular, and 
can be expressed as 



/37 



(56) 



, (57) 



where the symmetric matrices Fp^ and G^^ are given by 
the expressions 



1 



A 



(-) 



A, 



(-) 



(58) 



G 



A 



(+) 



A. 



(+) 



Pi 



ujp + uj^i + 2h{s + 1)1 crw^ + I (jujp + I 



(59) 



If uji is a real root, and uj2 



is a complex conjugate V. THE AMPLITUDE EQUATION 



pair, the equation for the phase shift can be written as 
follows 



B.B 



2»3 



aoji + I auj2 + I 



+ c.c. 



(52) 



Remarkably enough, in both cases the slowly varying 
phase of the basic harmonic s is a real number. Solv- 
ing now Eq. (1501) for the first-order correction to the 



The final step of our perturbation analysis consists in 
deriving an amplitude equation for the slowly varying en- 
velope functions Bj3, introduced in the zero-order approx- 
imation (|4T|) . For that purpose, let us write the second- 
order equation following from the expression p7p 



{adt + il){idt - hs f - ^{adt + *0 + ^ 
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a{idt — hs)^ + 2i{adt + il){idt — hs) — 



to the nonhnear wave- number shift in nonhnear optics^—. 
Note also that the second term on the right-hand-side of 
Eq. (|5T|) describes purely cross-modulation effects. In 
addition, the four-wave mixing effect responsible for en- 
ergy exchange between the three eigenmodes has been 
neglected. 



il) 



(idt 



(60) 



After substituting into the right-hand-side of the above 
equation the already found solutions in previous orders 
and performing simple but rather cumbersome algebraic 
manipulations, one can isolate the secular terms. These 
can be canceled by equating them to the second term on 
the left-hand-side of Eq. (15(11) . The result represents the 
amplitude equation governing the evolution of the slowly 
varying amplitudes Bp on a time scale of the order 1 / e^- 
times longer as compared to the fast oscillations with 
frequencies specified by the characteristic equation 

In the case, where all three roots of the characteristic 
equation (|^ are real, we obtain 



ihVf. 



7 = 1 



where 



(61) 



Vp = a{ujp + hsf + 2{aujp + l){ujp ^ hs) - (62) 



and 



aujp + I 
h{auj^ + 1) 



Tp^ = Fp^A\-^ - Gp^A\+^ 



[2{ujp + hs){Fp^ + Gp^) + h[Fp^ - Gp^)], 



(63) 



A^-^A^+^-A^^A^-^ 



A^+)^-A^-)^ 
aujp + I 



4:{aujp + I) / ujp + hs ujp 



+ 



(64) 



Similar amplitude equations can be obtained when one 
of the roots of the characteristic equation is real, while 
the other two comprise a complex conjugate pair. In this 
case, however, the exponential growth of the Fourier har- 
monics determined from zero and first order is dominat- 
ing and the nonlinear wave interaction yields relatively 
small contribution. 

The terms proportional to F^^ describe the effects of 
phase self- modulation, whereas those proportional to F^-^ 
(with /? ^ 7) are responsible for phase cross modula- 
tion. Both result in an intensity-dependent nonlinear fre- 
quency shift of the eigenfrequencies ujp, an effect similar 



VI. NONLINEAR STANDING WAVES 

Let us assume a trivial time dependence of the wave- 
function of the form 



(65) 



where ^(x) is a real valued function of the position vari- 
able. According to Eq. ([^0]) the single-wave amplitude a 
has a stationary solution, which can be expressed as 



a = i($2(2;)e-"). 



(66) 



Substituting the above expressions into the Schrodinger 
equation (fT8|) . we obtain 



-I- A$ + ^[e"($2(^)e-*«) + e-"($2(^)e*«)] $ = 0, 

(67) 

where 



A 



2lu 



M = 



(68) 



This latter equation describes the standing wave solu- 
tions of the time-dependent equation (jlSp . which are the 
states with definite energy, instead of a probability dis- 
tribution of different energies as discussed in the previous 
three sections. 

Taking into account the periodicity of the problem, it is 
natural to represent the solution of Eq. ([57)) as a Fourier 
series 



(69) 



similar to the representation given by Eq. (I?!]) . Thus for 
the Fourier amplitudes, we obtain 



\ n n / 

Although it is in principle possible to solve the nonlinear 
equations (j70p for a finite number of Fourier amplitudes, 
it is instructive to consider the perturbative solution dis- 
cussed below. 

The zero-order solution is obvious 



= AS„ 



(71) 
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where s is a particular Fourier harmonic number and A is 
an arbitrary complex constant. In first order, we obtain 

+ ^^'^>[Z^{A*^i%+A^i'>). (72) 

From the equation for ^i^-* it follows that Ai = and 
^i^'' = Ai, where Ai is a new constant, which can even- 
tually be set to zero. As before (see previous sections 
for comparison), the only nonzero Fourier harmonics ap- 
pearing in first order are the (s ± l)-st ones. These are 
linearly coupled via the equations 

(2. - 1 + f,\A\') + ^,AHi'_ll = 0, (73) 



/=-l,O=l,/!=l.lxl0 ' 



1 
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FIG. 1. Regions of complex s < 22 and real s > 22 roots 
of the characteristic equation (|42|l . The four- wave mixing 
(FWM) phase mismatch parameter 5a; = 2cj2 — oji — tJa, where 
iui < aj2 < is also shown. 



- fiA*Hi'}, + {2s + l- f,\Af) $ = 0. (74) 

The above system of linear homogeneous equations has 
a non trivial solution, provided its determinant is zero. 
This gives 



2^ 



(75) 



Clearly, the parameter must be negative, which is 
equivalent to the condition I = —1. For the (s ± l)-st 
Fourier amplitudes, we obtain 



€1 = 



(1) _2^iA^ 



{2s -If 



(76) 



where again B is an arbitrary complex constant. The 
second-order perturbation equation can be expressed as 

(^2_^2)^(2)_^^^(0) 



M<°lr( 



A*^ 



(2) 
s-1 



(2)* 
s+l 



(77) 



The nonlinear frequency shift A2 and the additional 
(s ± 2)-nd Fourier amplitudes appearing in second order 
can be found to be 



MB? 
(2.-1)2 



(78) 



(2) (2) ^^l^B*^ 

^+2 2(2s- l)(s + 1)' ^^-2 (2s-l)3(s-l)- 

(79) 

Note that the arbitrary to this end complex constants 
Ai and B can be used to satisfy the normalization condi- 
tion of the wave function. It is also important to mention 
at this point the following. The system ([7D|) comprises 
a set of n nonlinear equations for a finite number n of 
Fourier harmonics with A being kept as a free parame- 
ter. The latter can be determined from the normalization 
condition of the wave function. This usually leads to a 
nonlinear dispersion relation of n-th order, the roots of 
which determine the stability properties of the standing 
wave solution. 



VII. NUMERICAL SIMULATIONS 

In this section some of the predictions of the pertur- 
bation theory are benchmarked to the exact system dy- 
namics obtained by solving Eqs. (|18l) - (I^Hl) numerically 
by means of the standard split-step Fourier method—. 
Depending on the values of the relevant parameters the 
characteristic equation (|42l) can either posses three real 
roots or one real and two complex-conjugate roots. 

This is illustrated in Fig. [1] where the regions of 
three real and one real and two complex conjugate roots 
are marked for a particular choice of parameter val- 
ues (cr, I and h). In addition, in the region of real 
roots the value of the phase mismatch parameter 5ui = 
2a;2— wi— responsible for the four- wave mixing (FWM) 
proces&iS, is also shown. Within the limits of validity of 
the assumption for weak nonlinearity no evidence of en- 
ergy exchange between the three frequency components 
of either A^°_\ or A^^^ was found, which justifies neglect- 
ing the FWM in the derivation of Eq. (pT|) . This may be 
explained with the possible low value of the coefficient in 
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/=-l, a=l, /i=l.lxlO"^ j=10, e=10"' /=-l, a=\, /i=l.lxlO"^ j=10, e=10"' 




FIG. 2. Evolution of tlie modulus of tlie spectrum of tlie wave FIG. 4. Comparison between tlie time-domain numerical and 

function |vE'(s, t)\ with time. Initial conditions in the unstable analytical [Eq. (g!])] solutions for A''^^^ in the unstable region, 
region of parameter space have been chosen and interpolation 
has been used along the s-axis. 



fc-l,a=l,/7=l.lxlO -,j=48, £=5x10 




FIG. 3. The same as in Fig. [2]but the input is a stable modal 
triplet of the form (|80l) . 



front of the FWM term in the slowly- varying envelope ap- 
proximation (SVEA) equation (|6ip . A similar situation 
exists in the classical case of the single-wave modeiiS. 
As equations (|45|) and (j46|) show the amplitudes of the 

two weak sideband components and A^^^j^ either 

oscillate between their minimum and maximum values 
(three real roots), or alternatively, increase exponentially 
due to energy transfer from the powerful "pump" com- 
ponent \E'i°^ (two complex-conjugate roots of the char- 
acteristic equation). Thus modes excited in the region 
of purely real roots can be regarded as stable, whereas 
modes in the region of two complex-conjugate roots - as 
unstable. This observation is illustrated in Figures [5] and 

El 

Equations (IT51) - ^U\i have been solved with initial con- 
ditions in a triplet form 



*(a;,i = 0) 



1 



/27r 



Jsx ^ ^^^i{s~l)x ^i(s+l)x^ . (80) 



As Fig. [2] shows exciting an unstable modal triplet (with 
s < 22) leads to a gradual energy transfer to practically 
all neighboring modes and the excitation soon occupies 
the entire spectral range simulated numerically. In con- 
trast, in the stable region energy is initially transferred 
to neighboring modes but then returns back to the mode 
that has been excited originally and the process repeats 
itself periodically in time. 

In Figure |4] the analytical expression for the sideband 



component amplitude 



A 



(0) 



obtained by means of a 



perturbative approach is compared to the exact numer- 
ical result. The choice of parameter values corresponds 
to the unstable mode region. As can be seen excellent 
agreement between the two results exists at earlier times, 
whereas at later times they deviate from each other. This 
is because the perturbative solution (|45l) is valid pro- 



vided that 



A 



(0) 



(0) 



(weak probe interacting with 



a strong pump), which is indeed the case initially. How- 



ever, smce 



A 



(0) 



increases exponentially, the weak probe 

assumption is violated for sufficiently late times, which 
explains the observed difference between the analytical 
and the numerical result at later times. 

A similar comparison has been performed in the sta- 
ble region and the results are presented in both frequency 
(Fig. [S]) and time (Fig. ^ domains. A fast Fourier trans- 
form has been performed on A^'i\ = "^{s — l,t) and the 
result is presented by the full curve on Fig. [5j The three 
peaks correspond to the three real roots of Eq. (|42|) . In 
order to increase the frequency resolution a slow Fourier 
transform has then been performed in the vicinity of the 
peaks (dashed lines in Fig. [SJ. Finally, the amplitudes of 

the three components of A^'^^^^ have been computed from 
Eq. (|^5|) and plotted (the full circles in Figure [S]) . As 
can be seen the numerical and the analytical results are 
in excellent agreement. 



The same procedure has been applied to A 



(1) 



2,t) and the results are shown in Figures [7] and [SI Note, 
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FIG. 5. Numerically computed Fourier spectrum modulus 
|AfJi| (full and dashed curves) and (full circles) spectral am- 
plitudes obtained from Eq. (|45p in the stable region. Full 
curve: spectrum computed with the fast Fourier transform 
(FFT); dashed curves: spectrum computed with the "slow" 
Fourier transform (SFT). 
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FIG. 7. Frequency-domain solutions for jA^^jj- The curves 

represent numerically computed spectra (by FFT and SFT) 
and the full circles are the spectral amplitudes obtained ana- 
lytically from Eq. (f56)) . 
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FIG. 8. Time domain analytical and numerical solutions for 



FIG. 6. Time-domain numerical and analytical solutions for 
A'^I^j^ in the stable region. 



that according to Eq. (f56|) the spectrum of A^^j contains 
a total of six discrete lines corresponding to combination 
frequencies of the type ujfj + ujj in full agreement with 
what Fig. [7] shows. 

Figure [9] shows the evolution of the modulus of the 
spectrum of the wave function |\E'(m,i)| obtained from 
the numerical solution of Eq. (|18p with the amplitude 
of the single- wave mode a given by Eq. (|66p. The ini- 
tial condition 'ii{m,t — 0) is determined by Eqs. (1751) . 
(175)1 and ([75)1 . As can be seen, a stable solitary- wave- like 
structure is formed. 



o=l, /^=l.lxlO"^ 5=256, 5=10"^ 



100 



,240 220 

FIG. 9. Solitary-wave-like behavior of |vl/(m,t)|. The initial 
condition is given by Eqs. ([75]), ((76]) and ((79)l . 
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VIII. CONCLUDING REMARKS 

The present paper can be regarded as an extension 
of the classical single-wave model^ii^ to the case, where 
quantum-mechanical effects begin to play comparable 
role. The quantum single-wave model has been derived 
starting from the von Neumann-Maxwell equations and 
closely following the basic steps outlined by del-Castillo- 
NegreteS. It is important to note that the basic equations 
in Schrodinger representation figuring as a starting point 
for our analysis can be obtained formally by canonical 
quantization of the classical single- wave model. 

The linear stability of the quantum single- wave model 
has been studied in Sections IIIII and IIVI It has been 
shown that irrespective of the character of the spectrum 
of eigenfrequencies (real or complex) , once a certain har- 
monic is generated, it excites its closest neighbors due to 
the self-consistent closest neighbor interaction law, which 
is intrinsic for the system. In the stable region (of pa- 
rameter space, where all roots of the characteristic equa- 
tion are real) periodic in time patterns are formed, which 
have been observed to be marginally stable. In the unsta- 
ble region (one real and two complex conjugate roots of 
the characteristic equation) a numerically detected chain- 
type excitation of closest neighbors resembles a chaotic 
phenomenon. The freshly generated harmonics excite 
new ones, thus spreading gradually in harmonic number 
space until all the Fourier spectrum becomes populated. 

Based on the method of multiple scales, a system of 
coupled nonlinear equations for the slowly varying am- 
plitudes of interacting eigenmodes has been derived. As 
a matter of fact, the nonlinear standing-wave solutions 
of the single-wave model can be regarded as solutions 
of a nonlocal nonlinear Schrodinger equation.. A class 
of solutions has been found, which behave as a stable 
solitary-wave pattern. 



Appendix A: Linear Stability Properties of the Quantum 
Single-Wave Model 

The von Neumann equation coupled with the equation 
for the single field mode has an evident equilibrium solu- 
tion in the form / = /o(p), a = oq = 0. To analyze its lin- 
ear stability, we represent the Wigner quasi-probability 
distribution function as 



fix,p;t) ^ fo{p) + ^ fn{p]t)t 



(Al) 



and substitute it into Eq. (IT2]) . Retaining linear terms 
only, we obtain 



dtfn + inpfn = 



iaS„. 



±1 



/o ( P - 2 



{adt+il)a^i J dpfiip;t), (A3) 

where a is a deviation from the equilibrium value aq = 0. 
We further assume that the linear modes evolve according 
to the law 



/„(p; t) = F„(p)e^'^*, a{t) = ae"*^*, 
and readily obtain 



(A4) 



Fnip) 



aS. 



n,±l 



h{p - c) 



fo[p 



Mp 



2 J " V 2 
Finally, the dispersion relation can be expressed as 

dp/o(p) 



V{c) = ac-l~ J 



{p-cy - fi^/A 



0. 



(A5) 



(A6) 



Let us take the equilibrium Wigner quasi-probability dis- 
tribution in the form 



fo{p) = -^Sip-hs), 
Zn 



(A7) 



which is consistent with the exact solution (f26)) of the 
Schrodinger equation (fT8|) . The dispersion relation (jA6P 
transforms then to the form 



V{c) ^ac-l 



1 



1 



2n {hs-c)~tf/A 



= 0. 



(A8) 



(A2) 



Identifying now c with — w, the latter expression coincides 
with the characteristic equation (H^ obtained in Section 
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